Obtaining the size distribution of fault gouges with polydisperse bearings 
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We generalize the recent study of random space-filling bearings to a more realistic situation, where the spacing 
offset varies randomly during the space-filling procedure, and show that it reproduces well the size-distributions 
observed in recent studies of real fault gouges. In particular, we show that the fractal dimensions of random 
polydisperse bearings sweep predominantly the low range of values in the spectrum of fractal dimensions ob- 
served along real faults, which strengthen the evidence that polydisperse bearings may explain the occurrence 
of seismic gaps in nature. In addition, the influence of different distributions for the offset is studied and we find 
that the uniform distribution is the best choice for reproducing the size-distribution of fault gouges. 

PACS numbers: 46.55.+d, 45.70.-n, 61.43.Bn 



I. INTRODUCTION 

In the early nineties the question of the possibility to tile 
an arbitrarily large strip of space-filling roller bearings with- 
out friction nor slipping was addressed flfl, motivated by the 
study of real systems such as seismic gaps. Seismic gaps are 
regions along a fault zone where earthquakes do not take place 
and therefore they could be explained by sheared plates on a 
space-filling bearing @|. Fault zones are typically self-similar 
and the mechanical origin of the power-law in the particle size 
distribution was associated to the particle's fracture probabil- 
ity which has been proposed to be controlled by the relative 
size of its nearest neighbors iHQl- More recently, a geophys- 
ical model ISj] explained the different values of the fractal di- 
mension, ranging from fi/ = 2.6 to d/ ~ 3, by taking into 
account the fault gouge strain. While such model explains 
the dynamical origin of different power-laws, there is still the 
question if the space-filling bearing scenario is able to repro- 
duce such empirical results in a simple and systematic way. 
By reproducing with space-filling bearings the same particle 
size distributions observed in fault zones one can strengthen 
the hypothesis that the existence of seismic gaps in fault zones 
may be related to the emergence of particular geometrical ar- 
rangements of their composing rocks, due to local fragmenta- 
tion during tectonic motion. Figure [T] illustrates two random 
space-filling bearings in two and three dimensions. 

Pioneering studies with space-filling bearings were done 
using deterministic procedures in two ilf| and three di- 
mensions and also using random algorithms 10]. However, 
up to now specific initial configurations and fixed parame- 
ter values were addressed. In this paper we present a general 
algorithm to construct realistic space-filling bearings that al- 
lows to reproduce the range of fractal dimensions observed in 
fault gouges IH 3 [H Ht]. The model is parameterized by a 
unique parameter that controls the strength of fragmentation 
and takes into account ensemble averages. Despite the wider 
freedom in the parameters and initial configurations, the sys- 
tems presents robust results in what concerns the fractal di- 
mension. In particular, we will show that by varying the range 
of admissible values of the control parameter one finds fractal 
dimensions observed in fault gouges. 




FIG. 1 : (Color Online) Illustration of a random space-filling bearing 
in (a) three and in (b) two dimensions. Discs and spheres of the same 
color do not touch each other. The random space-filling bearing starts 
with a large disc or sphere, maximizing polydispersity (see text). 



We start in Section |ll]by describing in some detail the pro- 
cedure to generate random space-filling bearings, introducing 
a parameter that accounts for fragmentation at the local scale. 
In Sec. imiandllVlwe describe the results for two and three di- 
mensions respectively, with special emphasis on the size dis- 
tribution and the fractal dimension. Discussions and conclu- 
sions are given in Sec.lVl 
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n. THE RANDOM SPACE-FILLING OF PARTICLES 

In this Section we will start by revisiting previous proce- 
dures f7*| for constructing random space-filling packings and 
bearings and then introduce the necessary ingredients to ob- 
tain a fully random space-filling bearing. 

Random bearings in two and three dimensions are con- 
structed in the following way. First, one starts by randomly 
distributing a small number Nq of discs or spheres within a 
given range of sizes, without touching each other Second, 
one fills the empty spaces in the system by introducing itera- 
tively the biggest possible disc or sphere in the neighborhood 
of some empty region. Third, one resizes some disc or sphere 
in order for the packing to be bi-chromatic (bearing condi- 
tion), i.e. only two colors are needed to color all discs in such 
a way that no discs of the same color touch each other This 
guarantees the bearing condition: particles are able to role on 
each other without friction or slipping. Figures|2^ and|2j5 give 
illustrative examples of such random bearings in two dimen- 
sions. 

The filling procedure is done by choosing randomly a void 
within the inter-disc free space and then fitting the biggest disc 
in it, i.e. fit the disc that touches the three nearest discs in the 
neighborhood, as illustrated in Fig. |2^. For the three dimen- 
sional case one considers spheres touching the four nearest 
neighbors. 

The coloring procedure is done by attributing a proper color 
to the introduced disc. In the case that the three neighboring 
discs have the same color one attributes the other color to the 
new disc. Otherwise, one chooses only one of the neighbors 
to be in contact with the new disc, and the new disc shrinks 
to a size with radius r — aro (0 < a < 1), where vq is the 
radius before shrinking and gets a different color as the disk it 
touches. Figurellt illustrates this coloring procedure. 

Parameter a is our control parameter For constant a = 1 
one obtains the particular case where bearing cannot be guar- 
anted. In this case, frustrated contacts emerge when parti- 
cles are forced to rotate f^, which would eventually lead to 
the fragmentation of the discs into smaller ones. Recently 
a new method to implement realistic grain fracture in three- 
dimensional simulations of granular shear was proposed ||^, 
based on breakable bonds between particles within a medium. 
We keep the model simple, by using instead the reduction fac- 
tor a that mimics the effect of fragmentation: by shrinking a 
particle originally with frustrated contacts we mimic its frag- 
mentation into smaller particles that will fill the empty space 
left after the fragmentation. 

The algorithm described above was previously iFtI IToll used 
in two and three dimensions by fixing a given initial configu- 
ration with A^o discs having sizes within a given fixed range 
and also by using a fixed reduction factor a during the filling 
and coloring stages. The density of the packing was stud- 
ied as a function of the number N of existing spheres as 
well as the cumulative particle size distribution. It was found 
that the cumulative distribution obeys a power-law, namely 
N{r) = n{q)dq ^ r^'^f , where d/ is the fractal dimen- 
sion of the bearing U2R . 

Next, we introduce the additional points to strengthen pre- 



vious findings and improve the algorithm described above. 

First, there is the statistical significance of the results, and 
their sensitiveness to initial configurations, i.e. to the initial 
range [r* — S/2,r* + S/2] of sizes. This initial configura- 
tion may influence the polydispersity of the system and con- 
sequently the attained distribution after filling the entire sys- 
tem. As we will see, large number iVo of initial discs typically 
influence how the density increases during the space-filling 
procedure. 

To take this point into account we enable the construction 
procedure to start with an initial configuration having a single 
arbitrarily large disc (or sphere) and perform ensemble aver- 
ages on a significant number of initial configurations. Con- 
cerning the single initial large disc, one should notice that it 
does not suffice one single disc to introduce a second one, be- 
cause in two-dimensions each inserted disc needs to have at 
least three neighbors (four neighbours for three-dimensions). 
However, as illustrated in Fig.|2];, such starting disc can be in- 
troduced into the system by previously distributing a few very 
small 'seed-discs' in the system and then following the algo- 
rithm described above. The number of such seeds is small, 
and therefore they do not affect significantly the cumulative 
size distribution. Their role is that the average distance be- 
tween them is eventually of the order of the system size, en- 
abling the introduction of a first disc with the size of the order 
of the system size. Of course, depending on the number and 
distribution of the seeds, the first initial discs may have also a 




FIG. 2: (Color Online) Sketch of the construction of a random space- 
filling bearing, (a) A new disc (blue) is randomly inserted in the 
system and is shifted and enlarged to the maximal accessible size 
(dashed circle) without overlapping neighboring discs. Then, (b) it 
is reduced by a factor < a < 1, keeping a single contact point with 
one of the neighbors and assuming the opposite color, (c) To start the 
space-filling with a large disc (brown) one needs to place previously 
a few small 'seeds' (green) in the system (grey) and then proceed as 
in (a) and (b) [see text]. 
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FIG. 3: Density and size distribution in two-dimensional random space-filling packings, (a) Density p of the packing as a function of the total 
number of discs for different values of a = 0.2, 0.6 and 1.0, together with (b) the distribution of the radius r of the discs. In both cases, 
one starts from a fixed initial configuration of Nq — 40 discs with radius in the range [0.087?, 0.14i?] with R = 1 being the radius of the 
system (see Fig.[TJ)). The minimal radius of the initial set of discs is indicated as T^n 0.08 (see text). Fitting of the power-law range in (b), 
the distribution N{r) — hr^'^f yields the fractal dimension d/ as a function of a plotted in the inset. In (c) one plots the density as a function 
of iV for the initial sets of iVo = 40 discs in [0.05i?, O.IOJ?], [0.15J?, 0.207?] and [0.20i?, 0.257?], i.e. ranges with the same interval width but 
centered around different values, namely around r* = 0.075, 0.175 and 0.225 respectively and fixing a = 0.6. In the inset of (c) the density 
p{N) is plotted for initial ranges [r* — |, r* + |] having the same center r* = 0.075 but different widths 5 = 0.2, 0.4, 0.6 and 0.8. In (d) we 
plot distribution N{r-) of the spheres as a function of r for the same conditions as in (c). In all cases A'' = IC discs. 



size within the initial range of sizes. In this way, one general- 
izes the previous procedure 101 and maximizes the admissible 
polydispersity. 

Second, we also introduced a criterion to increase compu- 
tational efficiency of our algorithm. The neighborhood were 
the neighboring discs (or spheres) are searched for must be 
chosen conveniently. We propose to choose a size that de- 
creases with the increase of the density p, since the denser the 
packing the smaller are the empty spaces to put new discs. 
Therefore, the radius r„ of the neighborhood of a given ran- 
dom point introduced in the system at iteration n, is updated 

Tn = ^^^{rsys " r^ax), wherc rsys and r^ax are the 
radii of the system and of the biggest disc or sphere in it, re- 
spectively. 

Third, we also consider the control parameter a to vary ran- 
domly within a tunable range of values. In particular, we argue 
that though a tentative value of constant a could be obtained 
by analyzing samples of gouges in real situations, one expects 
that a certain range of admissible values for a is the most re- 
alistic assumption. Indeed, we show that the typical range of 



fractal dimensions observed in real fault gouges is in this way 
reproduced. 



III. THE TWO-DIMENSIONAL CASE 

We start this section by addressing the two-dimensional 
random space-filling bearing and systematically reviewing the 
behavior of the packing for different, but fixed, values of a 
and study the effect of fixed initial size ranges (no maximal 
admissible polydispersity). Figure |3] shows for this case the 
density p{N) and cumulative size distribution N{r) of two- 
dimensional random space-filling packings. 

Figure shows the density p as a function of the number 
N of discs for a = 1.0 (packing) and also for a = 0.2 and 
0.6 (bearings) separately, starting from TVo = 40 discs with 
radius in the range [r* — |,r* + |] having r* = 0.11 and 
S = 0.06. As expected, the convergence p 1 as increases 
is faster for larger values of a. We consider one fixed initial 
configuration with iVo initial discs, and therefore the different 
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FIG. 4: Averaging (a) the density p and (b) the distribution size N{r) over 100 initial configurations starting from large discs (r ^ R/2) and 
with a varying in a range [0.5 — Aa/2, 0.5 + Aa/2] with Aa — 0.2, 0.4, 0.6 and 0.8. The inset of (b) shows that the fractal dimension in 



N{r) 



is almost independent of Aa yielding dj ~ 1.54 which is, within the numerical errors of real fault gouges (see text). 



curves coincide for N < Nq. 

In Fig. ^ we plot the distribution of the radius r of the 



discs, where r„ 



is the minimal radius of the initial 



set of discs, and the deviation from the power law for r > ?■„ 
is due to the initial configuration. Below this value r„i, the 
size distribution obeys a power-law N — br~'^f , where d/ is 
the fractal dimension of the packing, plotted in the inset as a 
function of a (symbols). As one sees from the inset, the fractal 
dimension typically takes values in the range 1.2 < df < 1.4, 
differently from the values found in two-dimensional cuts of 
fault gouges 1.6 ± 0.1). There is a maximum of df for 
a = 0.5 that can be explained from the definition of a in the 
algorithm described above. For an a < 0.5, to each new disc 
introduced there is a remaining free space characterized by 
a' > 0.5 such that a + a' ^ 1, and similarly for a > 0.5. 

Both Figs.|3^ and|3}' consider the same initial configuration. 
To study the influence of the initial configurations, we plot in 
Fig. [St the density p{N), fixing a = 0.6, similar to previous 
works ifioll, and using different size ranges for the initial sets 
of A^o = 40 discs, namely in [OMR, O.IOR], [0.15R, 0.20R] 
and [0.20i?, 0.25i?], i.e. ranges with the same width S = 0.06 
but centered around different values, namely around r* = 
0.075, 0.175 and 0.225 respectively. 

Since different initial configurations are now used, the den- 
sity is no longer the same below A^o as in Fig. [3^. Further, one 
observes that the density converges to one for increasing the 
value of r*. In the inset the density p{N) is plotted by fixing 
r* = 0.075 and starting with initial configuration having dif- 



ferent widths, namely 6 = 0.2, 0.4, 0.6 and 0.8. In these cases 
the density gives always similar dependencies on A^. There- 
fore, the average size r* of the initial configuration is the im- 
portant parameter to tune the density of the packing. Its width 
can be varied without changing significantly the results. 

In Fig.|3}l we plot the accumulative size distribution N{r) 
of the discs for the same conditions as in Fig. The value 
of the exponent remains almost constant, df ^ 1.35. In other 
words, the fractal dimension is not very sensitive to the initial 
configuration, and the parameter on which the fractal dimen- 
sion depends more strongly must be indeed a. 

Since a is also the parameter controlling the fragmenta- 
tion of discs with frustrated contacts (see above), we will now 
study it more deeply. When a is able to vary randomly, the 
fragmentation of the largest disc in the free holes can be re- 
garded as a random process by its own. We next consider 
a to be each time randomly selected from a fixed interval 
[a* - Aa/2, a* + Aa/2]. We will show that when enabling 
a to take different values for each particle shrinking, one ob- 
tains fractal dimensions similar to the ones observed in fault 
gougesHiH. 

To this end, we put everything together, namely a varying 
in the middle range of admissible values, a large initial disc 
and an ensemble average over a significant number of initial 
configurations. The results for density and size distribution 
are shown in Fig.H] where one considers three initial seeds in 
the range ta/ = 2rm ^ i?/ 1000, with i? the size of the system 
and a varies randomly in the range [0.5 — Aq:/2, 0.5 + Aa/2]. 
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FIG. 5: The fractal dimension as a function of the exponent /3 when 
the value of a is chosen according to a power-law P{a) ~ in 
a range a £ [0.1, 0.9]. Although within the error bars, the fractal 
dimension is somewhat lower compared to Fig. |4j5 where the distri- 
bution of chosen a values is uniform (see text). 



Averages over a sample of 100 initial configurations. Since 
initialization, filling and coloring procedure are now all ran- 
dom, we call these systems fully random space-filling bear- 
ings. 

Figure 2^ shows the density as a function of the number TV 
of discs for Aa = 0.2,0.4,0.6 and 0.8. One sees an abrupt 
transition above = 3 (initial seeds), due to the introduction 
of the first large disc. For all the four cases the dependence 
of p on the range of a-values is similar, with the convergence 
towards p = 1 being slightly slower for narrower ranges, be- 
cause they hinder the occurrence of large discs. 

Figure I4J3 shows the size distribution N{r) for each of the 
four ranges. All the distributions almost coincide, as shown in 
the inset where the fractal dimension dj taken from N{r) ^ 
r~'^f is almost constant (df ^ 1.54). This value is larger 
than the one obtained when a is kept constant (see Fig. O. 
Notice that the value of the fractal dimension for Aa = 0, 
though corresponding to the case of constant a = 0.5, is dif- 
ferent from the one plotted in the inset of Fig. [3}', since the 
constructing procedure of the bearing is slightly different (see 
SecHD. 

Since the above value is obtained from a significantly larger 
sample of initial configurations and all the parameters a and 
position of the discs are randomly selected, we will consider 
this value df = 1.54 as the characteristic exponent of the 
size distribution for fully random two-dimensional space fill- 
ing bearings. The average characteristic value df = 1.54 ob- 
tained lies in the range of values measured of the fractal di- 
mension measured in real fault gouges (i? = 1.6±0.1)ili, 
as indicated with a dashed line and shadow region in the inset 
ofFig.Hh. 

Till now, the value of a was considered to vary uniformly 
within a certain range of values. If the probability distribution 
for choosing a values is Gaussian similar results are obtained. 



were the standard deviation <t of the distribution plays a simi- 
lar role as the width Aa used in Fig.|4] 

However, if we take values within a range, say a G 
[0.1, 0.9], and chose them according to a power-law distribu- 
tion P{a) ~ a^^ the fractal dimension changes significantly, 
as shown in Fig. |5] Even for small values of the exponent (3, 
e.g. (3 = 0.5, the fractal dimension decreases when compared 
to the value obtained for the uniform distribution {[3 = 0) and 
remains approximately constant at ^ 1.43. 

Notice that /? = in Fig. |5] corresponds to Aa = 0.8 
in Fig. |4}3. If the power-law distribution selects values in a 
range with a different width, a similar decrease of the frac- 
tal dimension is observed when comparing with the uniform 
distribution case. Therefore, one can conclude that a reason- 
able choice for constructing space-filling bearings with fractal 
dimension similar to the one observed in fault gouges is by 
taking a random value of a uniformly distributed in a certain 
range around 0.5. 



IV. THE THREE-DIMENSIONAL CASE 

As described above in Sec. [Ill a three-dimensional version 
of fully random space-filling bearings is obtained in a similar 
way as for discs with the single difference that the introduction 
of new spheres takes into account four nearest neighbors. In 
this Section we address the case of three-dimensional space- 
filling bearings as a more realistic approach to fault gouges, 
and study how well the two-dimensional model approximates 
three-dimensional systems of spheres. 

Recently ||^, it was found that grain fracture simulations 
produce a comminuted granular material similar to the one 
observed in real fault gouges. From those simulations, it 
followed that comminution rate and survival of large grains 
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FIG. 6: The fractal dimension as a function of a for three- 
dimensional random space-filling bearings when a is kept constant. 
A similar procedure as the one illustrated in Fig.|2]is used, with new 
spheres being introduced touching the four nearest neighbors (see 
text). The dashed line indicates one typical value df ~ 2.58 found 
in some real fault gouges iQ]. 
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FIG. 7: Three-dimensional space-filling bearings for varying a G [0.5 — Aa/2, 0.5 + Aa/2]. (a) The density p as a function of the number 
N of spheres for Aa = 0, 0.2, 0.6 and 0.8 and (b) the corresponding size distribution N{r) with the fractal dimension d/ in inset. In all cases 
N = 10*. 



is sensitive to applied normal stress, with a fractal dimen- 
sion of the resultant grain size distributions in the range 
df G [2.3 ± 0.3, 2.9 ± 0.5], that agrees with the observations 
of three-dimensional samples of real gouges where typically 
df ~ 2.58. 

In three dimensional space-filling bearings with a constant 
value of a the fractal dimension lies above the observed val- 
ues in fault gouges. In Fig.|6]we plot typical values of d/ as a 
function of a. The fractal dimension of such bearings is typ- 
ically larger than df > 2.58 (dashed line) with values within 
the range df e [2.60 ± 0.10, 2.74 ± 0.15] and, similarly to 
the two-dimensional case, the maximum of df is reached for 
a ~ 0.5. 

As summarized above in the introduction, it was recently 
found istl that fractal dimensions ~ 2.6 are observed for low- 
strain gouges. In regions subject to larger shear strain the frac- 
tal dimension is significantly larger, < 3. Therefore, the par- 
ticle size or mass dimensions were proposed as a way to dis- 
tinguish between regions with different strain strengths Js]]. 
From Fig. |6] one sees that a similar range of values for the 
fractal dimension is also found for space-filling bearings. 

Furthermore, the explanation relating the fractal dimension 
of fault zones and their strain strength assumes that fragmen- 
tation is controlled by nearest neighboring particle contact 
and that a particle is most likely to split into smaller parti- 
cles with a particle of similar size, yielding a larger fractal 
dimension ^ 3. In the case of our construction procedure 
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FIG. 8: (a) The size distribution of two-dimensional cuts of the three- 
dimensional space-filling bearings addressed in Fig. [7] and (b) the 
corresponding fractal dimension. 



for space-filling bearings this would correspond to the case of 
a ^ 0.5. Indeed, from Fig.|6]one observes that the maximum 
of the fractal dimension is reached for such a values yielding 
df = 2.74 ±0.15. 

When varying a randomly in a range around 0.5 and study 
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the dependence of the space-filUng bearing on the width Aa 
of the range [0.5-AQ;/2,0.5+Aa/2].InFig.|7}i one sees that 
the density increases faster for larger Aa, similarly to what 
was shown in Fig. |4^. As for the fractal dimension. Fig. [TJ? 
shows that it decreases slightly when compared with the case 
of constant a (Aa = 0). Therefore, increasing the width Aa 
of the range of admissible values for a one is able to reduce 
the fractal dimension of the bearing. 

Similarly to the situation of measures taken in fault gouges, 
the two-dimensional cross section of such three-dimensional 
space-filling bearings should have a fractal dimension within 
the range of the observed empirical values in real fault gouges 
(d/ = 1.6±0.lil|l). By averaging several different two di- 
mensional cross sections of the 3D bearings we plot in Fig.[8t 
the size distribution of a typical two-dimensional cross sec- 
tion for the different values of Aa. In Fig.[8}5 we observe that 
only for very wide ranges of a values it is possible to obtain a 
fractal dimension similar to the one observed on fault zones. 



V. DISCUSSION AND CONCLUSIONS 

In this work we studied the size-distribution of random 
space-filling bearings with large polydispersity, showing that 
it reproduces well the size-distribution found in fault gouges. 
Focusing on the dependence of the bearings fractal dimension 
on the spacing offset, we have shown that the fractal dimen- 
sions of such bearings sweep the low range of values observed 
in real faults. Since recently it has been reported that the frac- 
tal dimension varies in space along fault gougesisll, our find- 
ings enables us to conjecture that the occurrence of seismic 
gaps, where earthquakes are absent and therefore behave sim- 



ilarly to roller bearings, may occur in regions where the fractal 
dimension lies in the low range of admissible values, namely 
df G [2.5,2.75]. 

To compute an accurate value for the exponent characteriz- 
ing random bearings, we introduced a general algorithm that 
allows a to vary randomly in a wide range of admissible val- 
ues, typically < a < 1 and start the space-filling procedure 
from one unique large disc (or sphere), maximizing the range 
of admissible sizes in the bearing. 

With such model we were able to show that bearings have 
a fractal dimension with values within the range of values in 
real fault. Since it is known jSi] that along a specific fault 
gouge the fractal dimension varies typically between ~ 2 and 
~ 3, our results support the hypothesis that seismic gaps, oc- 
curring only in certain particular locations of the fault, could 
be explained by this simple geometrical model. 

Further, we also interpret the control parameter a for the 
bearing property as a measure of the fragmentation strength, 
and introduce simple criteria to improve the computational ef- 
ficiency of previous space-filling packing algorithms. 

To improve further our findings we should also take the ef- 
fect of gravity into account. Moreover, concerning the space- 
filling bearings by themselves, other questions raise, namely 
their contact network correlations, which should help to un- 
derstand the range of observed fractal dimensions. These and 
other points will be addressed elsewhere. 
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